Transient detection and flux measurement stand at the base of time domain astronomy. Due to changing seeing conditions, the detection of new sources on complex constant background is a highly non-trivial problem. Existing solutions for these problems experience large and significant problems, as hordes of image subtraction artifacts swamp real transients in numbers. The existing remedy for this problem uses highly complex solutions using machine learning algorithms to narrow the sea of candidates. Last, human scanners sift through the remaining candidates to extract the real transient signals from the stack of candidates.
Detection of previously unknown transient sources is at the base of many fields of astronomy. Examples for such are the searches for Supernovae and microlensing events. The problem of image subtraction has proven to be hard to tackle, with many suggested solutions e.g., [Phillips & Davis 1995, Astronomical Data Analysis Software and Systems IV, 77, 297; Alard, & Lupton, 1998, ApJ, 503, 325; Bramich, 2008, MNRAS, 386, L77; Gal-Yam, Maoz, Guhathakurta, & Filippenko, 2008, ApJ, 680, 550]. The reason for the difficulty in subtracting images is that the point spread function (PSF) of images taken from the ground is varying. This change in the observations' PSF can generate many subtraction artifacts in the naive subtraction of the images.
The solution used by most surveys is based on the algorithm of Alard & Lupton for the image subtraction itself. After that, a detection process occurs, such as applying a matched filter and searching for new sources. This solution is far from being perfect and many so called “subtraction artifacts” remain after its operation.
These subtraction artifacts that are produced by popular image subtraction methods such as Alard & Lupton and Bramich 2008] are caused by the need to solve a complex inversion problem. This inversion problem can be regarded as a regularization effort on the partial deconvolution done by Phillips & Davis. Apart of being slow, this inversion problem is in itself an effective deconvolution, as the numerical instability of the deconvolution process cannot be swept under the rug. In cases where the PSF of the new image is substantially different from the reference PSF, the above algorithms can either cause ringing artifacts to appear in the subtraction image due to the effective division in Fourier plane, or can fail to shape the PSF of reference image to the form of the new image PSF, leaving many subtraction artifacts of a naive subtraction with non-matching PSF's appear.
These artifacts prevent the automatic detection of transients by causing many false positive signals. The current state of the art solution to this problem is to train a machine learning algorithm, e.g., [Bloom, Richards, Nugent, et al. 2012, PASP, 124, 1175; Goldstein, D'Andrea, Fischer et al. 2015, arXiv:1504.02936; Wright, Smartt, Smith, et al. 2015, MNRAS] to filter most of the artifacts and reduce the number of false positives to the minimum. Later, human scanners sift through all remaining candidate detections and decide which is real and which is not.
This elaborated process can undermine the successful operation of a transient detection survey in many ways. The first is that employing many human scanners can be very cumbersome and expensive. Many current surveys are spending lots of manpower on candidate sifting, examples for such surveys are PTF [Law, Kulkarni, Dekany, et al. 2009, PASP, 121, 1395] and PANSTARRS [Kaiser, Aussel, Burke, et al. 2002, Proc. SPIE, 4836, 154]. This use of human scanners is un-scalable, and is unfeasible for future surveys like LSST [Ivezic, Tyson, Abel, et al. 2008, arXiv: 0805.2366].
The second is that it incurs time delay on transient detection, as it is very hard for human scanners to detect transients in real time. This can damage greatly many science cases in which it is of utmost importance to make rapid follow-up observations of new transients. Moreover, at least some machine learning algorithms throw away real obvious transients. Finally, the human scanning step makes it very difficult to estimate the completeness transient surveys as human scanners cannot be persistent in their decision making process. An example for such a science case is the characterization of the circumstellar environments of supernovae progenitors, which is possible only with rapid follow-up in timescales of hours Gal-Yam et al.
Another problem is that even human scanners can be unsure if a transient is real or not, and many surveys adapt the methodology of accepting only candidates that are persistent in two or more consecutive observations. Needless to say, this methodology trades survey speed with the credibility of candidates.
Last, this scanning methodology prevents transient detection at the faintest limit, as it makes the statistics of detection hard to quantify, and relies on visual inspection of human scanners which from obvious reasons cannot detect transients at the statistical detection limit.
Brief Overview of Existing Methods for Image Subtraction
The solutions suggested to image subtraction divide naturally to two variants. The first, and more popular, variant could in general be referred to as regularized partial deconvolution. Solutions in this family are of Phillips & Davis; Alard & Lupton; and Bramich 2008. Gal-Yam et al. suggested a second variant, which is called herein “symmetric matched filtering”.
Denoting the new image by N and its PSF by PN, the reference image by R and its PSF by PR, the first approach attempts to find a convolution kernel k such that:N−k⊗R≅0  (1)
The first solution suggested for finding such a kernel k was given by Phillips & Davis. They simply suggested to perform a deconvolution solution in Fourier space:
                              k          ^                =                                                            P                N                            ^                                                      P                R                            ^                                ≈                                    N              ^                                      R              ^                                                          (        2        )            where {circumflex over ( )} represents Fourier transform. However, this solution is numerically unstable as the deconvolution operation can (and many times does) involve division by small numbers.
Alard & Lupton suggested a practical way to solve the numerical instability problem. Representing k as a set of basis function, and noting that Eq. 1 is linear, they suggested solving for k using linear least squares. Alard & Lupton suggested using a set of basic functions, which are linear combinations of Gaussians multiplied with low degree polynomials. Later on, Bramich 2008 suggested solving for the values of a pixelized kernel. Both the Alard & Lupton method and the Bramich 2008 method can be viewed as regularization of the deconvolution method of Phillips & Davis, i.e. restricting the solutions for the kernel k to some set of “logical” solutions.
Even though the numerical stability of these algorithms is better than that of Eq. 2, they still have several problems: First, the division by zero problem is still there and it can become especially pronounced when the new image has a narrower PSF (actually, it suffices for it to be narrower in some axis) compared with the reference image (interestingly, these methods are a-symmetric to the new and reference image, while the problem is symmetric to the exchange of R and N). Second, there is no statistical justification for any of these methods, there is no rigorous proof they are not losing information, and it is unclear what further image processing should be applied. For example, does one need to apply another matched filter on the subtracted image? If so, which filter should one use? Third, there is no clear prescription on how to set a detection threshold for new source detection and how to decide if a detected source is real or an artifact, or to quantify the probability of it being a false positive.
The symmetric matched filtering solution suggested by Gal-Yam et al. is to match filter the new image with the PSF of the reference image, vice versa, and subtract. i.e.:SGal-Yam08=Pr⊗N−Pn⊗R  (3)
This solution is always numerically stable, and leaves no subtraction artifacts. The problem with this solution is that again, there is no statistical justification apart from common sense. This solution was generally disregarded in the community mostly due to the concern that it loses sensitivity. This claim can also be derived using the common sense in which if you apply it to two images with equal PSF's, the subtraction image has an effective PSF which is twice as wide as the original, leading to a loss in sensitivity.
It should be mentioned that the problem of subtracting two images, and minimizing the resulting difference image in the least square sense has an infinite number of solutions (also in [Yuan & Akerlof 2008, 677, 808]). For example, the linear equation:Kr⊗R−Kn⊗N≅0,  (4)where Kr and Kn are arbitrary kernels, has an infinite number of solutions. This is because for any Kr, one can find Kn that satisfies Eq. 4 in the least squares sense. This simple analysis shows that all the mentioned subtraction methods are focused on making the PSF of the two images identical, with very little attention to the maximization of the signal-to-noise ratio (S/N) of a transient source that appears in one of the images. In a sense, these methods do not solve the transient detection problem, but a different problem, which is how to make two images as similar as possible, using convolution.
Accordingly there is a long felt need for a method that can be numerically stable, maximally sensitive and fast to compute.